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Abstract —In this paper, we analytically study the performance 
of multiple-input multiple-output underwater wireless optical 
communication (MIMO UWOC) systems with on-off keying 
(OOK) modulation. To mitigate turbulence-induced fading, which 
is amongst the major degrading effects of underwater channels 
on the propagating optical signal, we use spatial diversity 
over UWOC links. Furthermore, the effects of absorption and 
scattering are considered in our analysis. We analytically obtain 
the exact and an upper bound bit error rate (BER) expressions 
for both optimal and equal gain combining. In order to more 
effectively calculate the system BER, we apply Gauss-Hermite 
quadrature formula as well as approximation to the sum of 
lognormal random variables. We also apply photon-counting 
method to evaluate the system BER in the presence of shot 
noise. Our numerical results indicate an excellent match between 
the exact and upper bound BER curves. Also a good match 
between the analytical results and numerical simulations confirms 
the accuracy of our derived expressions. Moreover, our results 
show that spatial diversity can considerably improve the system 
performance, especially for channels with higher turbulence, e.g., 
a 3 x 1 MISO transmission in a 25 m coastal water link with 
log-amplitude variance of 0.16 can introduce 8 dB performance 
improvement at the BER of 10 -9 . 

Index Terms —Underwater wireless optical communications, 
MIMO, spatial diversity, lognormal fading channel, photon¬ 
counting approach, saddle-point approximation, optimal combin¬ 
ing, equal gain combiner. 

I. Introduction 

U NDERWATER wireless optical communications 
(UWOC) has recently been introduced to meet a 
number of demands in various underwater applications due 
to its scalability, reliability, and flexibility. As opposed to 
its traditional counterpart, namely acoustic communication, 
optical transmission has higher bandwidth, better security, 
and lower time latency. This enables UWOC as a powerful 
alternative for the requirements of high-speed and large-data 
underwater communications such as imaging, real-time video 
transmission, high-throughput sensor networks, etc. 0-0. 
However, despite all the above advantages, due to the severe 
degrading effects of the UWOC channel, namely absorption, 
scattering, and turbulence 0, currently UWOC is only 
appropriate for short-range communications (typically less 
than 100 m) with realistic average transmitted powers. 
Removing this impediment, in order to make the usage 
of UWOC more widespread, necessitates comprehensive 
research on the nature of these degrading effects and also 
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on the intelligent transmission/reception methods and system 
designs. 

In past few years, various studies have been carried out, 
theoretically and experimentally, to characterize the absorption 
and scattering effects of different water types. Mathematical 
modeling of a UWOC channel and its performance evaluation 
using radiative transfer theory have been presented in [8j. A 
novel non-line-of-sight UWOC network has been proposed 
in 0, based on the back-reflection of the propagating op¬ 
tical signal at the ocean-air interface. In 1101, a hybrid op¬ 
tical/acoustic communication system has been developed. The 
beneficial application of error correction codes in improving 
the reliability and robustness of UWOC systems has been 
experimented in ED- Based 011 experimental results in 0 
and 1121 for absorption and scattering of different water types, 
the beam spread function and the channel fading-free impulse 
response were characterized in d and 0, respectively. 
Furthermore, a cellular underwater wireless network based on 
optical code division multiple access (OCDMA) technique has 
been proposed in 0. while the potential applications and 
challenges of such an underwater network and the beneficial 
application of serial relaying on its users’ performance have 
very recently been investigated in d and d respectively. 

Unlike acoustic links where multipath reflection induces 
fading on acoustic signals, in UWOC systems optical turbu¬ 
lence is the major cause of fading on the propagating optical 
signal through the turbulent seawater d- Optical turbulence 
occurs as a result of random variations of refractive index. 
These random variations in underwater medium mainly result 
from fluctuations in temperature and salinity, whereas in at¬ 
mosphere they result from inhomogeneities in temperature and 
pressure changes d’d- Recently, some useful results have 
been reported in the literature to characterize the underwater 
fading statistics. A precise power spectrum has been derived 
in d for the fluctuations of turbulent seawater refractive 
index. Based on this power spectrum, the Rytov method has 
been applied in 0, p0| to evaluate the scintillation index of 
optical plane and spherical waves propagating in an under¬ 
water turbulent medium. Also Tang et al. d have shown 
that temporal correlation of irradiance may be introduced by 
a moving medium and they investigated temporal statistics of 
irradiance in moving ocean with weak turbulence. Furthermore 
in d, the on-axis scintillation index of a focused Gaussian 
beam has been formulated in weak oceanic turbulence, and by 
considering a lognormal distribution for intensity fluctuations, 
the average BER has been evaluated. Moreover, the average 
BER and scintillation index of multiple-input single-output 
(MISO) UWOC links have very recently been analyzed in |[22j 
and |23J, respectively, when degrading effects of inter-symbol 
interference (ISI) are neglected. Additionally, experimental 
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Fig. 1. Block diagram of the MIMO UWOC system with OOK modulation. 

studies on the fading statistics of UWOC channels in the 
presence of air bubbles have quite recently been carried out 
in (24). 

Despite all of the valuable research have been done on 
various aspects of UWOC, a comprehensive study on the 
performance of UWOC systems that takes all degrading effects 
of the channel into account is missing in the literature. While 
some work only considered absorption and scattering effects 
and used numerical simulations, without analytical calcula¬ 
tions, to estimate the BER of point-to-point UWOC systems 
[1], some others considered turbulence effects, but neglected 
scattering and ISI effects, and used the common conditional 
BER expressions in free-space optics (FSO) pT) , (25) . The 
research in this paper is inspired by the need to comprehen¬ 
sively evaluate the BER of UWOC systems with respect to 
all impairing effects of the channel. Moreover, we employ 
multiple-input multiple-output (MIMO) transmission to inves¬ 
tigate how MIMO technique can mitigate turbulence effects 
and extend the viable communication range. The technique of 
spatial diversity, i.e., exploiting multiple transmitter/receiver 
apertures (see Fig. [T]), not only compensates for fading effects, 
but can also effectively decrease the possibility of temporary 
blockage of the optical beam by obstruction (e.g., fish). An¬ 
other advantage of spatial diversity is in reducing the transmit 
power density by dividing the total transmitted power by the 
number of transmitters. In other words, depending on the 
wavelength used there exist some limitations on the maximum 
allowable safe transmitted power. Hence, by employing spatial 
diversity one can increase the total transmitted power by the 
number of transmitters and therefore support longer distances, 
while maintaining the safe transmit power density 09- 

In this paper, we apply maximum-likelihood (ML) detection 
to analytically obtain the BER expressions for MIMO UWOC 
systems with repetition coding across the transmitter^] when 
either an optimal combiner (OC) or an equal gain combiner 
(EGC) is used. In order to take the absorption and scattering 
effects into account, we obtain the channel impulse response 
using Monte Carlo (MC) simulations similar to (T), (26). 
Moreover, to characterize the fading effects, we multiply the 
above impulse response by a fading coefficient modeled as a 
lognormal random variable (RV) for weak oceanic turbulence 
(2T), (23). Closed-form solutions for the exact and upper 

'Note that because of the positive nature of incoherent wireless optical 
communications, a non-destructive addition of different transmitters’ intensity 
signals appears at the receiver; therefore, full channel diversity can be achieved 
by using simple repetition coding, i.e., sending the same on-off keying (OOK) 
symbol from all of the transmitters for each bit interval. 


bound BER calculation are provided using the Gauss-Hermit 
quadrature formula as well as approximation to the sum of 
lognormal RVs. We also apply the photon-counting approach 
to evaluate the BER of various configurations in the presence 
of shot noise. In our numerical results we consider both 
spatially independent and spatially correlated links. 

The rest of the paper is organized as follows. In Section 
II, we review some necessary theories in the context of 
our proposed MIMO UWOC system and the channel under 
consideration in this paper. In Section III, we analytically 
obtain both the exact and upper bound BER expressions, 
when either OC or EGC is used at the receiver side. We 
also apply Gauss-Hermite quadrature formula to effectively 
calculate multi-dimensional integrals with multi-dimensional 
finite series. In order to evaluate the system BER using photon¬ 
counting methods, the same steps as Section III are followed in 
Section IV. Section V presents the numerical results for various 
system configurations and parameters considering lognormal 
distribution for fading statistics of the UWOC channel. And 
Section VI concludes the paper. 

II. Channel and System Model 

In this section, we present the channel model that includes 
all three impairing effects of the medium, followed by the 
assumptions and system model that we have introduced in 
this paper. 

A. Absorption and Scattering of UWOC Channels 

Propagation of optical beam in underwater medium induces 
interactions between each photon and seawater particles either 
in the form of absorption or scattering. Absorption is an irre¬ 
meable process where photons interact with water molecules 
and other particles and lose their energy thermally. On the 
other hand, in the scattering process each photon’s transmit 
direction alters, which also can cause energy loss since fewer 
photons will be captured by the receiver aperture. Energy 
loss due to absorption and scattering can be characterized by 
absorption coefficient a( A) and scattering coefficient b( A), re¬ 
spectively, where A denotes the wavelength of the propagating 
light wave. Moreover, total effects of absorption and scattering 
on the energy loss can be described by extinction coefficient 
c (A) = a (A) + 6(A). These coefficients can vary with source 
wavelength A and water types (T|. It has been shown in (5), 
(27) that absorption and scattering have their minimum effects 
at the wavelength interval 400 nm < A < 530 nm. Hence, 
UWOC systems apply the blue/green region of the visible light 
spectrum to actualize data communication. 

In Q,(26) the channel impulse response has been simulated 
based on MC method with respect to the absorption and scat¬ 
tering effects. We also simulate the channel impulse response 
similar to Q, @ relying on MC approach. In this paper, 
the fading-free impulse response between the ith transmitter 
and the jth receiver is denoted by ho,ij(t). As it is elaborated 
in 0, when the source beam divergence angle and the link 
distance increase the channel introduces more ISI and loss 
on the received optical signal. On the other hand, increasing 
the receiver field of view (FOV) and aperture size increases 
the channel delay spread while decreasing its loss. It is worth 
mentioning that although this behavior may be unobservable in 
clear ocean links, as the water turbidity increases and multiple 
scattering dominates it becomes more apparent. 
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B. Fading Statistics of UWOC Channels 

In the previous subsection we described how the absorption 
and scattering effects are characterized. To take turbulence 
effects into account, we multiply ho ij (t) by a multiplicative 
fading coefficient a'X [18], [281—[32 , with lognormal distri¬ 
bution for weak oceanic turbulence |2T). (23). The assump¬ 
tion of lognormal distribution for UWOC turbulence-induced 
fading statistics is a reasonable assumption. In fact, UWOC 
systems typically employ large receiving apertures due to the 
negligibility of background noise under water which leads to 
an effectively reduced turbulence. Furthermore, although the 
attenuation due to pointing error and misalignment for FSO 
links can be modeled as a RV, for the sake of simplicity and 
owing to the highly-scattering nature of UWOC channels, the 
perfect alignment is also assumed in this paper. To model 
turbulence-induced fading let a = exp (X) be the fading 
amplitude of the channel with lognormal probability density 
function (PDF) as 


fa (a) = 


1 




i exp 


(In (a) -tixf 
2 (J 'x 


( 1 ) 


Therefore, the fading log-amplitude X has a Gaussian dis¬ 
tribution with mean yx and variance a 2 x . To ensure that 
fading neither amplifies nor attenuates the average power, we 
normalize fading amplitude such that E [a 2 ] = 1, which 
implies fix = [ 331. In order to thoroughly describe 

fading statistics of UWOC channels we should find out the 
dependency of log-amplitude variance to the ocean turbulence 
parameters. 

The scintillation index of a light wave is defined by 0. 

(D; 


a 2 (r, d 0 , A) 


(I 2 (r,d 0 ,\)) - (I (r, d 0 , A)) 2 
(/Mo, A)) 2 


( 2 ) 


in which I(r,dp,X) is the instantaneous intensity at a point 
with position vector (r, do) = (x,y,dp), where do is the 
propagation distance and (•) denotes the long-time average. 
Assuming weak turbulence, the scintillation index of plane 
and spherical waves can be obtained as 0- ® 


pi POO 

aj = 8n 2 kld 0 / / k $„( k ) 

Jo Jo 

dpk 2 


1 0 -to 
X { 1 — cos 


£( 1 - ( 1 - 0)0 




dndf, (3) 


in which 0 = 1 and 0 for plane and spherical waves, 
respectively, kp = 2tt/X and k denote the wave number 
and scalar spatial frequency, respectively; and <1>„ (n) is the 
power spectrum of turbulent fluctuations given by Eq. (8) of 
|7J. For lognormal turbulent channels the scintillation index 
relates to the log-amplitude variance as a 2 = exp(4<r^-) — 1 
p8 |. In [|7|, pO) the scintillation index of plane and spherical 
waves is numerically evaluated versus communication distance 
do for various values of scintillation parameters and it is 
observed that, depending on the values of these parameters, 
the strong turbulence (which corresponds to a 2 > 1 0 , @) 
can occur at distances as long as 100 m and as short as 10 m, 
which impressively differs from atmospheric channels where 
strong turbulence distances are on the order of kilometers 0- 
Therefore, mitigating such a strong turbulence demands more 
attention. 


C. MIMO UWOC System Model 

We consider a UWOC system where the information signal 
is transmitted by M transmitters, received by N apertures, and 
combined using OC/EGC. As it is depicted in Fig. [T] optical 
signal through propagation from the ith transmitter TXi to 
the jth receiver RXj experiences absorption, scattering, and 
turbulence, where in this paper absorption and scattering 
are modeled by fading-free impulse response of ho t ij(t ) (as 
discussed in Section II-A) and turbulence is characterized by 
a multiplicative fading coefficient a 2 ,-, which has a lognormal 
distribution for weak oceanic turbulence (as explained in 
Section H-B). Therefore, the TXi to RXj channel has the 
aggregated impulse response of hij(t) = a 2 jhp t ij(i). 

We assume intensity-modulation direct-detection (IM/DD) 
with OOK signaling where the “ON” state signal will be 
transmitted with the pulse shape P(f). In the special case 
of rectangular pulse, P(t) can be represented as P ( t ) = 
PII where P is the transmitted power per bit “1”, 

T b is the bit duration time, and 11(f) is a rectangular pulse 
with unit amplitude in the interval [—1/2,1/2]. Therefore, the 
transmitted signal can be expressed as; 

OO 

S(t)= Y hP(t-kT b ), (4) 

k= — oo 

where b k £ {0,1} is the fcth time slot transmitted bit. 
Hence, the received optical signal after propagating through 
the channel can be represented as; 

OO 

y (t) = S (t) * a 2 h 0 (t.) = Y b k a 2 T{t-kT b ) , (5) 

k— —oo 

in which hp(t) is the channel fading-free impulse response 
and r (f) = hp(t) * P(f), where * denotes the convolution 
operator. 

When M transmitters are used, in order to have a fair 
comparison with the single transmitter case, we assume that 
the total transmitted power for the “ON” state signal is yet 
P = p, where Pi is the transmitted power by r fX,. 

In this regard, we denote the transmitted pulse from TXi by 
Pi ( t ) = PiU anc | a j so a sequence of data signal 

from TXi is denoted by S t (t) = J2kL-oo bkPi (t - kT b ). 
This signal reaches RXj after propagation through the channel 
with impulse response a 2 jho t ij(t). Hence, the transferred 
optical signal from TXi to RXj will be as; 

OO 

Vi,j ( t) = Si (f) * a 2 jh 0 ,ij(t ) = Y hc^jTij (f - kT b ) , (6) 

k— —oo 


where Ti j(t) = Pi (f) * hp^j (f). At the jth receiver, transmit¬ 
ted signals from all transmitters are captured, each with its own 
channel impulse response. In other words, received optical 


signal at the jth receiver is as yj (f) = Y^iLiVij (t) 118 
(30). Moreover, in order to have a fair comparison between 
multiple receivers and single receiver cases, we assume that 
the sum of all receiving apertures’ areas in multiple receivers 
scheme is equal to the receiver aperture size of the single 
receiver scheme. 

We should emphasize that all of the expressions in the 
rest of this paper are based on Fj ; j(f), which is in terms 
of Pi (f) and hpjj (f). Therefore, our analytical expressions 
are applicable for any form of power allocation between 
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different transmit apertures, any pulse shape of transmitted 
signal Pi (t), and any channel model. However, our numerical 
results are based on equal power allocation for transmitters, 
i.e.. Pi = P/M , i = 1,..., M, rectangular pulse for OOK 
signaling, MC-based simulated channel impulse response, and 
lognormal distribution for fading statistics. Additionally, all 
of the derivations throughout the paper are for a general 
case from the transmitters and receivers structures point of 
view, i.e., all of the transmitters and receivers are located in 
arbitrary places. This generality is covered by considering a 
specific impulse response for any transmitter-to-receiver pair 
as hij(t ) = ctijh.o t ij(t). However, in our numerical results 
we assume that all of the transmitters and also all of the 
receivers are located with an equal separation distance in a 
line perpendicular to the transmission direction. 

III. BER Analysis 

In this section, we analytically derive the exact and an upper 
bound BER expressions for both single-input single-output 
(SISO) and MIMO schemes, when either OC or EGC is used. 
Various noise components, i.e., background light, dark current, 
thermal noise, and signal-dependent shot noise all affect the 
system performance. Since these components are additive and 
independent of each other, in this section we model them as 
an equivalent noise component with Gaussian distribution (34). 
Also as it is shown in [|35j, the signal-dependent shot noise has 
a negligible effect with respect to the other noise components. 
Hence, it is amongst the other assumptions of this section to 
consider the noise variance independent to the incoming signal 
power. Moreover, we assume symbol-by-symbol processing at 
the receiver side, which is suboptimal in the presence of ISI 
[361. In other words, the receiver integrates its output current 
over each Tj, seconds and then compares the result with an 
appropriate threshold to detect the received data bit. In this 
detection process, the availability of channel state information 
(CSI) is also assumed for threshold calculation. 


threshold, i.e., with T = a 2 "/ 8 ' 1 /2. Therefore, the conditional 
probabilities of error when bits “ 1 ” and “ 0 ” are transmitted 
can respectively be obtained as; 


pSISO 

■ f foe|l,a:,bfc 


Pr(4iso < T\b 0 = 1 ) 


Q 


7 (s) /2 + El ; 


a T b 



( 8 ) 


pSISO 

■ f be|0,a',6fc 


Pr(4iso > T\b 0 = 0) 


Q 


7 (s) /2-El ; 


a T b 



(9) 


where Q (x) = (1/\/2 tt) exp{—y 2 /2)dy is the Gaussian- 

Q function. Then the final BER can be obtained by averaging 
the conditional BER P^ a ° bk = ^+ 5 ^° over 
fading coefficient a and all 2 L possible data sequences for b^s 
as; 


p be so = 4 E f p hY /«(«)**• a°) 

1 b k 

The forms of Eqs. © and ([9]> suggest an upper bound on 
the system BER, from the ISI point of view. In other words, 
bk^t 0 = 0 maximizes while Q has its maximum value for 
bk^o = 1- Indeed, when data bit “0” is sent, the worst effect 
of ISI occurs when all of the surrounding bits are “1” (i.e., 
when bkjto = 1), and vice versa [36|. Regarding these special 
sequences, the upper bound on the BER of SISO UWOC 
system can be evaluated as; 


pSISO 

r k,UB 



Q 


/«V s) \ 

V 2 VT b ) 


+ 


Q 


(a 2 

7 (s) /2-E^L7 (/ ’ fc) 

)} 


( J T b 

J 


f a (a)da. 


( 11 ) 


A. SISO UWOC Link 

Based on Eq. ([5|, the photodetector’s 0th time slot integrated 
current in SISO scheme can be expressed as; 

r siso = &oa 2 7 (s) + a 2 51 b kl P ' k) + v Tb , (7) 

k=—L 

where 7 W = R T(t)dt, 7 ( 7 > fc ) = R T(t - kT b )dt = 
R f_j,T b 1 ^ Tb T (t)dt, and R = r/q/hf is the photodetector’s 
responsivity. Moreover, 77, q, h, /, and L are the photodetec¬ 
tor’s quantum efficiency, electron’s charge, Planck’s constant, 
optical source frequency, and channel memory, respectively. 
Physically, 7 ( J > fe #o) re f ers t 0 the ISI effect and j( I ’ k = 0 ') 
interprets the desired signal contribution, i.e., 7 ( / > fc= °l = 7 O). 
Furthermore, vr b is the receiver integrated noise component 
which has a Gaussian distribution with mean zero and variance 
Oj, [341. Note that based on the numerical results presented 
in fl7| , the channel correlation time is on the order of 10 7 ‘ to 
10 " 2 seconds which implies that thousands up to millions of 
consecutive bits have the same fading coefficient. Therefore, 
we have adopted the same fading coefficient for all of the 
consecutive bits in Eq. ID- 

Assuming the availability of CSI, the receiver compares its 
integrated current over each Tj seconds with an appropriate 


The averaging over fading coefficient in m and CD 
involves integrals of the form J 0 °° Q(Ca 2 )f a (a)da, where 
C is a constant, e.g., C = 7 ^) /2crT b in the first integral of 
CD- Such integrals can effectively be calculated using Gauss- 
Hermite quadrature formula [37, Eq. (25.4.46)] as; 


Q(Ca 2 )f a (a)da 

/ OO 

Q(Ce 2x ) 

-OO 


y/2* 


exp - 


{x - yxf 
2a 2 x 


dx 


9=1 


-/=^2,w q Q ( Cexp ( 2x q j2a 2 x +2p x ) ) , ( 12 ) 


in which U is the order of approximation, w q , q = 1,2,..., U, 
are weights of the C/th-order approximation and x q is the 7 th 
zero of the [/th-order Hermite polynomial, Hu{x) |l 8 |, |37). 


B. MIMO UWOC Link with OC 

Relying on (| 6 j, the integrated current of the jth receiver can 
be expressed as; 


M 


M 


Jbo) 

3 


= b ° +Y, a l Y. b 


(i,k) 

,3 


o) 

V T b > 


(13) 
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where 7 $ = R f^ b r itj (t)dt, j£ k) = Rf^ b Tij(t- 

kTb)dt = R f_kT 1 ' >Tb L,, j is the memory of the 

b . • (i) 

channel between the zth transmitter and yth receiver, and Vj-* 

is the jth receiver integrated noise which has a Gaussian 
distribution with mean zero and variance a4- . It is worth 

-L b 

mentioning that while EGC simply adds the output of each 
receiving branch (with equal gain) to construct the combined 
output, the OC applies maximum a posteriori (MAP) detection 
rule (which in the case of symbols with equal probability 
simplifies to ML rule) to obtain the optimal combining policy. 

Using the derived decision rule in Appendix A, we find 
the “ON” and “OFF” states conditional error probabilities as 
Eqs. ( [T4| and ( p~5j ), respectively, shown at the top of the next 
page. Assuming the maximum channel memory as L max = 
iriax{L-|-|. L 12 , Lmn}, the average BER of MIMO UWOC 
system can be obtained by averaging over fading coefficients 
vector a (through an (M x AO-dimensional integral) as well 
as averaging over all 2 imax possible sequences for />y, : s, i.e., 

pMIMO _ 1 / 1 r dMIMO , 

“be ,OC 9 _L max 2.^ U O U 6e,OC|l,<3,bfc + 

b k Ja 

p ™!oc\o,a,b k ]fM)da, C16) 


where f&(a) is the joint PDF of fading coefficients in <5. 

Furthermore, considering the transmitted data sequences as 
bk^L 0 = 1 f°r bo = 0 and b^o = 0 for b 0 = 1 , the upper 
bound on the BER of MIMO UWOC system can be obtained 
as Eq. shown at the top of the next page. Moreover, as 
it is shown in Appendix B, (M x A')-dimensional integrals 
in © and © can effectively be calculated by (M x AO¬ 
dimensional series using Gauss-Hermite quadrature formula. 

It is worth mentioning that for transmitter diversity (N = 1) 
the conditional BER expressions in ( [T4| and © simplify to; 



tS+(-i ) 6 o+l E 2 ^; 


(im 

1 


which can be 

(-i) b0+1 E ~ k L L 


k=-Li 


reformulated 


2 &T b 


(18) 


as 
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where G 


(bo) 


pMISO^ 

J be\b 0 ,di,bk 


2,1 


b® 


a il G iA 

]kX-L u 2b kli*r']/2<7T h - The weighted sum 
of RVs in ( p~ 8 | ) can be approximated by an equivalent 
RV, using moment matching method (38) . Therefore, 
the conditional BER of dTSli 


we can 


in 
using 

approximate 


pMISO^ 

1 be\b 0 ,oi,bk 


q{ 


G 


(bo)\ 


as 


M J ’ 


M is the equivalent RV 


where G (bo) 

resulted from the approximation to the weighted sum of RVs, 
i.e., G^ « Efci a ii G i^i' > ■ In the special case of lognormal 
fading, the equivalent lognormal RV, G = exp( 22 l b °l), 
has the log-amplitude mean and variance of; 


MzO 0 ) 


M 


= - 


a %o) = 7 ln 


j=l 


1 + 


2 


E^i (g£ 0) ) 2 (< 


Mi 


il - 1 


(sr^M s~i(bo 
1 ^ 2,1 


(bo)\ 2 


(19) 


( 20 ) 


respectively (33) . Hence, in the case of transmitter diversity 
the average BER can approximately be evaluated with a one¬ 
dimensional integral which can also be calculated using Gauss- 
Hermite quadrature formula as a one-dimensional series. 


C. MIMO UWOC Link with EGC 


When EGC is used, the integrated current of the receiver 
output, based on Eq. ©, can be expressed as; 

N M N M -1 

r Mmo= +EXXX! bk ^ k) + v n ] > 

j = 1 2 — 1 j = 1 2 = 1 k— — Lij 

( 21 ) 

where v is the integrated combined noise component which 
has a Gaussian distribution with mean zero and variance Na^ b 


Based on © and the availability of CSI, the receiver 
selects the threshold value as T = VE, EE, 7 E/ 2 . 
Pursuing similar steps as Section III-B results to \22\ for 
the conditional BER. As expected, ( |22) simplifies to ( f]~ 8 j ) for 
MISO scheme. Finally, the average BER can be evaluated 
similar to ©• Also the upper bound on the BER of MIMO 
UWOC system with EGC can be expressed as Eq. ( f23) , shown 
at the top of the next page. 

It is worth noting that the numerator of ( |22] > can be 
approximated as <E o) « EjLi EEi E-.jXr where 
the weight coefficients are defined as Z)E = if- + 
(_l)b 0 + 1 Y2k=~L i:j ‘Zbk'Yij ^■ Similar to ( p~9| > and ( |20) , statis¬ 
tics of the equivalent lognormal RV C, lb °K which is resulted 
from weighted sum of M x N RVs, can be obtained and 
then averaging over fading coefficients reduces to the one¬ 
dimensional integral of; 

£( h o) 


Pi 


MIMO 


be,EGC|b 0 ,b fc 


Q 


2 y/N, 


a T b 


W C (bo) R (bo) , 


(24) 


which can also effectively be calculated using Eq. ©• 


IV. BER Evaluation Using Photon-Counting 
Methods 


In this section, we derive the required expressions for the 
system BER using photon-counting approach. Moreover, in 
this section signal-dependent shot noise, dark current, and 
background light all are considered with Poisson distribution, 
while thermal noise is assumed to be Gaussian distributed 
(36]|. To evaluate the BER, we can apply either saddle-point 
approximation or Gaussian approximation which is simpler 
but negligibly less accurate than saddle-point approximation. 
Based on saddle-point approximation the system BER can be 
obtained as Pb e = ^ [< 7 + (/3) + < 7 _(/?)], in which q+(/3) and 
g_(/3) are probabilities of error when bits “0” and “1” are 
sent, respectively, i.e.. 


q+ (/3) = Pr (zi > j3\zero) 

q- (/3) = Pr (u < (3\one) f 
$b 0 (s) = In (s)] - 


_ exp [$ 0 (sq)] 

\/ 27 r$ o( s o) ’ 
o exp [$i(ai)] 

s/3 — In |s|, b 0 


= 0,1, 


(25) 


where u is the photoelectrons count at the receiver output and 
'P u (6 0 ) (s) is the receiver output moment generating function 
(MGF) when bit “ bo ” is sent. Also so is the positive and real 
root of < f>p(s), i.e., ^(so) = 0 and Si is the negative and real 
root of $ 1 (s), i.e., c E > 1 (si) = 0; and /3 is the receiver optimum 
threshold and will be chosen such that it minimizes the error 
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probability, i.e., dP be /d(3 = 0. As an another approach to 
evaluate the system BER, Gaussian approximation is very fast 
and computationally efficient, yet not as accurate as saddle- 
point approximation, but yields an acceptable estimate of the 
system error rate particularly for BER values smaller than 0.1 
(36) . Indeed, when the receiver output is as u = N + £, 
where AT is a Poisson distributed RV with mean for the 
transmitted bit “&o” and £ is a Gaussian distributed RV with 
mean zero and variance er 2 , Gaussian approximation which 
approximates AT as a Gaussian distributed RV with equal mean 
and variance results to the following equation for the system 
BER (36); 

/ m (!) - m (0) \ 

Pbe = Q ( , f1 . , 2 - . 2 ). (26) 

V v Tfi'P + a 2 + v iw°> + a 2 J 

In this section, the required expressions for both the exact and 
the upper bound BER evaluations using either saddle-point or 
Gaussian approximation are presented, when EGC is used at 
the receiver side. 

A. SISO Configuration 

Based on Eq. 0 , the photo-detected signal generated by 
the integrate-and-dump circuit of the SISO receiver can be 
expressed as; 

u siso = Vsiso + v th, (27) 


where v t h corresponds to the receiver integrated thermal 
noise and is a Gaussian distributed RV with mean zero and 
variance a 2 h = 2K b T r T b /(RLq 2 ). where K b , T r , and Rl are 
Boltzmann’s constant, the receiver equivalent temperature, and 
load resistance, respectively [391. Conditioned on {bk}^l_L 
and a, y^i^O * s a P°i ss °n distributed RV with mean TOgtsc, 
as; 

2 0 rTb 

TO siso = TT E bk / r ( i - kT b) dt+(n b +n d )T b , (28) 

h J k——L Jo 


in which n b and n d are mean count rates of Poisson distributed 
background radiation and dark current noise, respectively. 

As it is shown in Appendix C, conditioned on a the receiver 
output MGF can be expressed as; 


i£ s ) =ex P 

“sisco 1 * 


2_2 


s (t: 


<ri 

k=—L 


th ' TO siso + b 0 a 2 m^ 


1 + exp ( e s — 1)) 


(e s - 1 ) 
(29) 


in which mgjgg = (n b + rid)T b , = £ j f^ b r(t)dt, and 

= Ff fo Tb r (^ - kT »)dt = § bn 1)Tb T ^ dt More - 

over, assuming the transmitted data sequences as bk^o = 1 for 
b 0 = 0 and b^o = 0 f° r = 1, MGF of the receiver output 
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for evaluation of upper bound on the BER of SISO UWOC 
system can be obtained as; 


•a. * <•>-•»«■ 


m siso + b 0 a 2 m^ 


+ ^2 b o a 


2 _(/,&) 
m y ’ ’ 


k— — L 


(e s -l) , (30) 


integrated thermal noise of the jth receiver and y^ condi¬ 
tioned on {bk}k=-L-- an ^ i * s a P°i sson distributed 

RV with mean; 

M 0 


m 


M 0 «T b 

(M = 777E E a l b */ r i,i (t-m)dt+(n dij +n b>j )T b , 


i=l k=-L 


(31) 


computationally fast. It can easily be shown that conditioned 
on bk the receiver output signal is the sum of a Gaussian 
and a Poisson RVs; therefore, Gaussian approximation can be 
applied to evaluate the average BER conditioned on a and bk, 
i.e., Pb e \a,b k - The Gaussian distributed RV has mean zero and 
variance Na/ h . And the Poisson distributed RV has mean of; 


where bo = 1 — 60 . Inserting ( |29[ i and ([30]) in ( |25[ > results into 
the conditional BER, P be \ a , and the final BER can then be 
obtained by averaging over the fading coefficient a. 

B. MIMO Configuration with EGC 

In this scheme, each of N receiving apertures receives 
the sum of all transmitters signals. At the receiver side, 
each of these N received signals passes through its receiver 
photodetector and different types of noises are added to each 
output. Therefore, the photo-detected signal at the jth receiver 
generated by integrate-and-dump circuit can be expressed as 
1 = yj b °' > +v th ,j, where v t h,j is a Gaussian distributed RV 
with mean zero and variance cr/ h = a/ h corresponding to the 


„_(6o) __(6rf) 

m MIMO — m MIMO 


N M 

j=i *=1 


boa: 
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2 ■TYl''- 

W 1 L i,j 


J2 bka t 


2 (I,k) 

iMj 


k—-Li 


(34) 


Moreover, mean of the Poisson distributed RV for the evalu¬ 
ation of the upper bound on the BER of UWOC system can 
easily be obtained by assuming the transmitted data sequences 
as bk^o = 1 for b 0 = 0 and b k ^o = 0 for bo = 1 . 

Using ([26]) and the conditional BER, P be \a,b k , can 
easily be evaluated based on Gaussian approximation. More¬ 
over, to obtain P be \s.,b k based on saddle-point approxima¬ 
tion, the simplified form of saddle-point approximation [36, 
Eqs. (5.73)-(5.79)] can be applied to ( |34| i. Subsequently, 
Pbe\b k can be obtained through an (M x iV)-dimensional 
integration as P be \ bk = f & P be \a,b k fa(a)da which yet 
demands excessive computational time, especially for large 
number of links. Nevertheless, we can reformulate ([34]) as 


(bd) 


“ '"-MIMO 

EEi 


' N 

b k m { */\ 
‘MIMO 


■ i\ ^ M 

2^7=1 2 ^ i=i 


T (bo) a 2 
i,3 iJ 


where t} 


»>j 


ml 


Hence, we can approximate 

= m MIMO + ^ (bt>) ’ Where 


where n b> j and n,i.j are the mean count rates of Pois¬ 
son distributed background radiation and dark current noise 
of the jth receiver, respectively. As it is demonstrated 
in Appendix D, MGF of the receiver output in MIMO 
scheme conditioned on fading coefficients vector a = 
(an,ai 2 , ...,aMJv) can be expressed as Eq. ([32|p| in which 
to mimo = {n b + Nn d )T b , m[*] = ft f^ b T zJ (t)dt, and 

m { -f = fiS?T itj (t-m)dt = fiS2£ m riAt)dt. 

Furthermore, MGF of the receiver output for the evaluation of 
the upper bound on the BER of MIMO UWOC system can 
be expressed as Eq. shown at the top of the next page. 

We should emphasize that extracting the output MGFs for 
MISO and single-input multiple-output (SIMO) schemes is 
straightforward by respectively substituting N = 1 and M = 1 
in ( [32| and Note that for these cases = m^ 0 

and TOgj ^ 0 = Eventually, using saddle-point ap¬ 

proximation the conditional BER P be \a can be achieved by 
inserting ( |32| ) and ( [33} in ( |25] >. The final BER can then be 
evaluated by averaging over a as P be = P be \afa (a) da. 

With respect to the above complex expressions, using 
saddle-point approximation for BER evaluation may be dif¬ 
ficult and computationally time-consuming, since it needs to 
solve some complicated equations for which their complexity 
increases as ISI (or equivalently L max in ( [32] )) increases. But 
( |26| suggests that using Gaussian approximation is simple and 

-Note that each of the receivers introduces a Gaussian distributed thermal 
noise with mean zero and variance a/ h . = nfi and a Poisson distributed 
dark current with mean count of jT b = n,jT b . But mean of the Poisson 
distributed background noise is proportional to the receiver aperture size and 
we assumed that the sum of all receiving apertures is identical to the aperture 
size of MISO scheme, which implies that jj ' V | n b j = n b . 


TO ( b °) 

"‘MIMO 

borriij t z^fc=-i, 

El ns Wmimo 

^mimo * s tb e approximated version of ( f34] > and $( b °) ss 
Ef=i r E° a ij , be., the weighted sum of M x N RVs. 
In the special case of weak oceanic turbulence the equivalent 
lognormal RV, d( b °) = exp(2 z^), has the following log- 
amplitude mean and variance, respectively 
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(36) 


By means of the above approximation, P b e\b k can 
be evaluated through two-dimensional integral of 

Pbe\b k « C CPbelb^^K^^d^d^, 

where P be \ bk is the system BER conditioned on b k , 

i9(°) and and (when Gaussian approximation is used) 

can be obtained as; 





~(1) ~(0) 
TO MIMO m MIMO 


mimo+ No th + y/rhlmio+NatM 

(37) 


Note that in the MISO scheme all of the transmitters are 
pointed to a single receiver; therefore, all of the links have 
the same fading-free impulse response and channel mem¬ 
ory as /lo.Miso(i) an d £miso> respectively. Consequently, 
when all transmitters have identical transmitted power of 
P/M, all links have equal m| s J and mf^ as to^ so and 
to miso , respectively. Hence, we can rewrite ( [34] ) as Wmiso = 
to miso + (bom-^so + Ea=Ll M iso bk <£)^ where 
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TABLE I 

Some of the important parameters used for noise 

CHARACTERIZATION AND MC-BASED CHANNEL SIMULATION. 


Coefficient 

Value 

Quantum efficiency, 77 

0.8 

Optical filter bandwidth, A A 

10 nm 

Optical filter transmissivity, Tp 

0.8 

Equivalent temperature, T e 

290 K 

Load resistance, Rl 

100 n 

Dark current, I^ c 

1.226 x 10“ a A 

Receiver half angle FOV, #fov 

40 u 

MISO schemes aperture diameter, 

20 cm 

Source wavelength, A 

532 nm 

Water refractive index, n 

1.331 

Source full beam divergence angle, 

0.02 u 

Photon weight threshold at the receiver, w t h 

10 -b 

Separation distance between the transmitters and 
between the receiving apertures, Iq 

25 cm 


¥> (M) = i- e -> the sum of M lognormal 

RVs. As a result, P be | bfc for MISO scheme can be 
evaluated through one-dimensional integral of Pbe\b k ~ 
/q 30 Pbe\b k ,v( M '>f^ M ^)d { P < ' M ' > - Then, if the channel memory 
is Tmax bits, Pb e can be obtained by averaging as Pb e = 
2 zh^'£b k p be\b k - Note that as ISI increases, this averaging 
demands more computational time and evaluation of the upper 
bound BER becomes more advantageous. 

From ( [34] ), one can observe the destructive effect of ISI on 
the BER. In other words, experiencing more time spreading 
in ho.ij(t) or equivalently Tjj(f) increases XAo m iE an£ I 
decreases m\ s ^. Thereby, it causes an increase in and a 
decrease in which results into larger BERs. Constructive 
effect of spatial diversity is appeared as combining the fading 
coefficients of different links which can be approximated as 
a single lognormal RV with roughly a scaled log-amplitude 
variance by the number of links | [TF| . 

V. Numerical Results 

In this section, we present the numerical results for the 
BER performance of UWOC systems in various scenarios. 
We consider lognormal distribution for the channel fading 
statistics, equal power as P/M for all transmitters, the same 
fading statistics (log-amplitude variance) for all links and the 
same aperture area of A/N for all of the receivers, where 
A is the total aperture area. In simulating the turbulence- 
free impulse response by MC method, we consider coastal 
water and also turbid harbor water links with absorption and 
scattering coefficients of (a, b) = (0.179,0.219) m” 1 and 
(0.366,1.824) m -1 , respectively |6j. Other important parame¬ 
ters for MC simulations are listed in Table I. In addition, some 
of the important parameters for characterization of noises are 
addressed in this table and the other parameters are exactly 


the same as those mentioned in ®, HD- Based on these 
parameters, noise characteristics are as rib ~ 1.8094 x 10 8 
s -1 in 30 meters deep coastal water, rid ~ 76.625 x 10 8 s _1 , 
and <y/ h /Tb = 3.12 x 10 15 s -1 . Hence, background radiation 
has a negligible effect on the system performance. 

Moreover, in this section we obtain the BER values using 
numerical simulations to verify the accuracy of our derived 
expressions. To do so, we generate 10 7 Bernoulli RVs with 
PDF of Pr(&o) = ^3(bo) + ^$(bo ~ 1), where S(.) is Dirac 
delta function. Subsequently, 10' pulses of shape boPi(t) will 
be transmitted to the channel. Each of these 10' pulses will 
be convolved with the channel fading-free impulse response 
and multiplied by a lognormal RV to constmct boa/jPi(t) * 
ho t ij(t) = boafjTij{t). Additionally, a contribution of chan¬ 
nel ISI reaches the receiver as described in Eqs. ([7]) and 
At the receiver side, the integrated current will be added by 
a Gaussian distributed noise and the result will be compared 
with an appropriate threshold (as comprehensively studied in 
Section III) to detect the received signal. Finally, comparing 
the transmitted data sequence with the detected one determines 
the BER value in a specified amount of the transmitted power. 

In order to see how water turbidity affects the channel loss 
and temporal dispersion, and also to investigate the channel 
spatial beam spread, we use MC method to simulate the 
channel fading-free impulse responses for a 1 x 2 SIMO 
transmission in a 25 m coastal water link and also an 8 m 
turbid harbor water link. The transmitter is pointed to the 
first receiver and the link between the transmitter and the 
first receiver has the impulse response of h 0 ,n(t). The second 
receiver is located in 25 cm center-to-center distance from 
the first aperture and receives those photons that are reached 
with much more scattering; with impulse response of ho , 12 (f)- 
Each receiving aperture has a diameter of 20/v/2 cm. Fig. 2 
illustrates the simulation results. Comparing Figs. 2(a) and 
2(c) shows that while the direct link in a 25 m coastal water 
has a negligible scattering, an 8 m turbid harbor water link 
remarkably scatters propagating photons and induces much 
more delay spread on ho, 11 (f)- In other words, as the water 
turbidity increases both the channel loss and delay spread 
considerably raise. As a result, an 8 m turbid harbor water link 
has a poor channel condition even than a 25 m coastal water 
link. Moreover, comparing Figs. 2(a) and 2(b) demonstrates 
that the second receiver mainly receives those photons that 
have experienced more scattering than the received photons 
by the first aperture. Accordingly, the channel delay spread in 
ho, 12 (f) is by far more than h 0 , 11 (f)- 

Fig. 0 depicts the exact BER of a 25 m coastal water 
link with transmitter diversity and data transmission rate of 
Rb = 1 Gbps. This figure also indicates an excellent match 
between the results of analytical expressions and numerical 
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Fig. 2. Fading-free impulse responses for a 1 X 2 SIMO UWOC system in different water types, (a) /io,n(i) of a 25 m coastal water link; (b) /io, 12 (f) of 
a 25 m coastal water link; (c) ho, 11 (f) of an 8 m turbid harbor water link. 




Fig. 3. Exact BER of a 25 m coastal water link with SISO, 2x1 MISO 
and 3x1 MISO configurations, obtained using both analytical expressions 
and numerical simulations. B;, = 1 Gbps, ax = 0.1 and 0.4. 

simulations. Here, we assume that the fading of each link is 
independent from the others. As it is obvious, increasing the 
number of independent links provides significant performance 
improvement in the case of ax = 0.4, e.g., one can achieve 
approximately 6 dB and 9 dB performance improvement at the 
BER of 10 l2 , using two and three transmitters, respectively. 
But this benefit relatively vanishes in very weak fading condi¬ 
tions, e.g., ax = 0.1. This is reasonable, since in very weak 
turbulence conditions fading has a minuscule effect on the 
performance but scattering and absorption have yet substantial 
effects. Hence, in such scenarios multiple transmitters scheme, 
which combats with impairing effects of fading, does not 
provide a notable performance improvement. 

In Fig. [4] we assume the same parameters as in Fig. [3] 
and use our derived analytical expressions to evaluate the 
exact BER of 1 x 2 SIMO, 1x3 SIMO and 2 x 2 MIMO 
configurations with optimal/equal gain combiner. Comparison 
between the results shows that the performance of EGC is 
very close to the performance of OC receiver. Therefore, due 
to its lower complexity, receiver with EGC is more practically 
interesting. Furthermore, the good match between the analyti¬ 
cal results and numerical simulations confirms the accuracy 
of our derived analytical expressions for the system BER. 
Comparing the results of BER for SISO and SIMO schemes 
demonstrates that a 1 x N' SIMO scheme provides better 
performance than alxJV SIMO configuration (A/ 7 > N > 1) 
only at high signal-to-noise ratios (SNRs) or equivalently low 


Fig. 4. Exact BER of a 25 m coastal water link with 1x2 SIMO, 1x3 SIMO 
and 2x2 MIMO configurations and optimal/equal gain combiner, obtained 
using both analytical expressions and numerical simulations. B;, = 1 Gbps, 
ax =0.1 and 0.4. 

BERs, where fading has more impairing effect than absorption 
and 1ST This is reasonable, since each receiver in a 1 x N' 
SIMO scheme has N' /N times less aperture area than a lxA f 
SIMO configuration and also an IV'-receiver scheme imposes 
N'/N times more dark current and thermal noise. In low 
SNRs, absorption and scattering as well as noise have more 
dominant effects on the BER than fading; therefore, in low 
SNRs the 1 x A" SIMO scheme yields better performance 
than the 1 x N' SIMO structure. However, when the channel 
suffers from relatively notable turbulence, the 1 x N' SIMO 
structure which has more links can better mitigate fading and 
can compensate for the loss due to the smaller aperture size 
and excess noise and therefore can yield better performance 
at higher SNRs. Needless to say that a 2 x 2 MIMO structure 
has the same aperture size as a 1 x 2 SIMO structure and 
since benefits from more independent links can yield better 
performance, than a 1 x N' SIMO, in all ranges of SNR. One 
can expect that in a very weak turbulence scenario, such as 
ax = 0 . 1 , dividing the receiver aperture to extend the number 
of independent links can degrade the performance. 

In Fig. [5] we applied our derived analytical expressions to 
evaluate the exact and upper bound BERs of a 25 m coastal 
water link with Rb = 1 Gbps and ax = 0.4, using (M x AO¬ 
dimensional integrals. As it can be seen, the upper bound 
BER curves have good tightness with the exact BER curves. 
Therefore, the upper bound BER evaluation can be more 
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Fig. 5. Comparison between the exact and the upper bound BERs of a 25 m 
coastal water link with R 5 = 1 Gbps, ax =0.4, and various configurations. 
Also the upper bound BERs are calculated with (M x A^)-dimensional series, 
using Gauss-Hermite quadrature formula (GHQF). 



Fig. 6 . Comparing Gaussian and saddle-point approximations in evaluating 
the exact BER of a 25 m coastal water link with R 5 = 1 Gbps, ax =0.4, 
and various configurations. Also the results of photon-counting methods are 
compared with the results of our derived analytical expressions for the system 
exact BER. 

preferable since the exact BER calculation may need excessive 
time for averaging over bk s. Furthermore, the upper bound 
BERs of various configurations are calculated by (M x AO- 
dimensional series, using Gauss-Hermite quadrature formula 
(GHQF). The order of approximation U is assumed to be 
the same for all of the links, i.e., Uij = 30. The excellent 
match between the results of GHQF and numerical (M x N)- 
dimensional integrals demonstrates the usefulness of GHQF 
in effective calculation of the system BER. 

Fig. [6] compares Gaussian and saddle-point approximations 
in evaluating the exact BER of a 25 m coastal water link with 
Rb = 1 Gbps, ax = 0.4, and various configurations. It is 
observed that Gaussian approximation can provide relatively 
the same results as saddle-point approximation. Therefore, 
due to its simplicity and acceptable accuracy, Gaussian ap¬ 
proximation can be considered as a reliable photon-counting 
method for the system BER evaluation. Moreover, the results 
of our derived analytical expressions are compared with those 
of photon-counting methods. The good match between the 
results of analytical expressions and photon-counting methods 
further confirms the validity of our assumption in neglecting 


the signal-dependent shot noise in our analytical derivations. 

After confirming the accuracy of our derived analytical 
expressions, a comparison between the computation time of 
different methods would be interesting. Table II shows the 
computation time and BER values of different methods for 
various system configurations in a 25 m coastal water link 
with ax = 0.4, Rb = 1 Gbps, and average transmitted power 
per bit of 20 dBm. The values in this table are obtained through 
(M x AQ-dimensional integration over fading coefficients as 
well as averaging over all 2 imax sequences for bk s. In order 
to perform the integration, we have taken A samp equidistant 
samples from the interval onj £ [0.0001,5]. Here, we have 
chosen L max = 3, A r salnp = 30 for MIMO configuration, 
and A^ sam p = 100 for the other schemes. As it can be seen, 
saddle-point approximation has by order of magnitude larger 
computation time, since it involves solving some complex 
nonlinear equations to obtain the conditional BERs. On the 
other hand, our derived analytical expressions and Gaussian 
approximation have relatively the same calculation time, since 
both involve Gaussian-Q function which is a fast-calculable 
built-in function in Matlab. Moreover, all three methods yield 
similar results for the system average BER, confirming the 
accuracy of our derived analytical expressions and the validity 
of our assumption in Section III, i.e., the negligibility of 
signal-dependent shot noise. A more advantage of our derived 
analytical expressions becomes visible when we apply the 
approximation to the sum of RVs. In this case, as it is 
elaborated in Section III, in the case of transmitter diversity 
(N = 1) and also receiver diversity with EGC, the (M x N)- 
dimensional averaging integration can be evaluated through an 
equivalent one-dimensional integral using our derived analyti¬ 
cal expressions, while Gaussian approximation only in the case 
of transmitter diversity leads to a one-dimensional integral, as 
shown in Section IV. Accordingly, when approximation to the 
sum of RVs is applied, our derived analytical expressions can 
more quickly predict the system performance. Therefore, our 
analytical approach is more advantageous and convenient from 
the computation time and ease of mathematical manipulation 
points of view; and hence is recommended to be employed for 
analysis of MIMO UWOC systems. 

In Fig. [7] the BER performance of a 25 m coastal water link 
with Rb = 0.5 Gbps and ax = 0.3 is depicted for different 
configurations. Also the sum of independent lognormal RVs 
in ( |34| ) is approximated with a single lognormal RV and 
the BER is evaluated through the approximated one-or two- 
dimensional integrals. As it can be seen, relatively good match 
exists between the results of the approximated one-or two- 
dimensional integrals and the exact (M x A)-dimensional in¬ 
tegral of Pbe = f s Pbe\ 3 .f& ( a)da. However, the discrepancy 
increases when receiver diversity is used. Moreover, the BER 
performance of a SISO link with ax =0.1 is compared with 
the BER curve of a 9 x 1 MISO link with ax = 0.3, and 
approximately the same result is observed; therefore, spatial 
diversity manifests its effect as a reduction in the fading log- 
amplitude variance. 

In order to investigate how ISI affects the performance of 
different configurations, in Fig. 8 the upper bound on the BER 
of a 25 m coastal water link with ax = 0.4 and various 
configurations is illustrated for a typical data rate of Rb = 0.5 
Gbps and also an extremely large data rate of Rb = 50 
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TABLE II 

The computation time (T CO mp) and BER values of different methods for various system configurations. Average transmitted 

POWER PER BIT = 20 DBM, R b = 1 GBPS, ax = 0.4, AND do = 25 M COASTAL WATER. 


___Configuration 

Method ~ ——___________ 

SISO 

2x1 MISO 

1x2 SIMO 
with OC 

1x2 SIMO 
with EGC 

2x2 MIMO 
with OC 

2x2 MIMO 
with EGC 

Analytical 

Tcomp 

0.05750 s 

1.5257 s 

1.7164 s 

1.6162 s 

140.6670 s 

127.8600 s 

Pbe 

3.050 x 10 -4 

5.088x10-° 

3.916X10- 4 

6.408 x 10 -4 

2.188x 10 t> 

2.639X 10- tt 

Gaussian approx. 

Tcomp 

0.05717 s 

0.8511 s 


0.8668 s 


116.3134 s 

Pbe 

3.056 x 10 -4 

5.106x10-° 


6.372 xl0“ 4 


2.607x lO -5 

Saddle-point approx. 

Tcomp 

38.4672 s 

1969.026 s 


7799.913 s 


432614.65 s 

Pbe 

3.004x 10 -4 

5.139x10-° 


6.288 x 10 -4 


2.639X 10“ s 



Fig. 7. Upper bound and approximated upper bound on the BER of a 25 
m coastal water link with different configurations, obtained using Gaussian 
approximation. R b = 0.5 Gbps, ax = 0.3 and 0.1. 



Fig. 8. Effect of ISI on the performance of a 25 m coastal water link with 
different configurations, ax = 0.4, R b = 0.5 Gbps and 50 Gbps. 

Gbps. As it can be seen, increasing the transmission data rate 
considerably degrades the performance of all configurations, 
especially for multiple receivers schemes where a lower bound 
on the BER is observable. This behavior can thoroughly be 
inferred from Figs. 2(a) and 2(b). In other words, in MISO 
schemes where a single receiver is employed, all of the 
transmitters are pointed to a single receiver; hence, all of 
the links from each transmitter to the receiver have impulse 
response of /io, 11 (f), shown in Fig. 2(a). In fact, based on Fig. 
2 (a), in both of the used bit duration times, i.e., Tf, = 2 x 10” 9 
s and 2 x 10 -11 s, the direct link delay spread is not much 


more than the bit duration time to see a saturation behavior in 
the BER of MISO schemes. Nevertheless, as the bit duration 
time decreases (data rate increases) the destructive effect of 
ISI becomes more severe and performance of both SISO 
and MISO schemes degrades; however, MISO transmission 
in both of the aformentioned data rates can yet yield a better 
performance than SISO transmission. On the other hand, in 
SIMO and MIMO schemes where more than one receiver 
is used, each transmitter is pointed to one of the receivers; 
however, due to the multiple scattering of UWOC channels 
some of the transmitted photons by each transmitter may reach 
the other receiving apertures that are not necessarily pointed to 
that transmitter. In other words, in SIMO and MIMO schemes, 
there always exist photons that have reached the receiver 
plane through indirect path, i.e., with the indirect link impulse 
response, such as /io, 12 (f) shown in Fig. 2(b). As it has been 
shown in Figs. 2(a) and 2(b), for Rb = 0.5 Gbps delay spread 
of both the direct and indirect links is enough smaller than the 
bit duration time; hence, ISI does not considerably degrade 
the performance of SIMO and MIMO schemes with Rb = 0.5 
Gbps. Nevertheless, for Rb = 50 Gbps, although the direct 
link temporal spread is yet enough smaller than the bit duration 
time, the indirect link delay spread is much more than the bit 
duration time. Therefore, for the case of Rb = 50 Gbps the 
destructive contribution of indirect link extremely increases 
ISI for multiple receivers schemes, where even increasing the 
transmission power cannot compensate for this degradation. In 
other words, in very large data rates the destructive effect of 
indirect link saturates the BER of multiple receivers schemes 
to a lower bound. It should be noted that decreasing the 
separation distance between the receiving apertures alleviates 
this time spreading, but this approach increases the spatial 
correlation between different links’ signals, which degrades 
the system performance (as it is shown in Fig. 11). In addition, 
we can reduce the ISI effect by decreasing the receiver aperture 
size, but this scheme induces more loss on the received optical 
signal. 

Fig. 0 demonstrates the ability of spatial diversity in ex¬ 
tending the viable communication range of UWOC systems, 
as one of the important objectives of this paper. In this figure 
we want to show that although increasing the communication 
distance increases absorption (loss), scattering (ISI) and tur¬ 
bulence (fading), yet using spatial diversity one can achieve 
better performance than SISO links with smaller distances 
and therefore can increase the viable communication range. 
To show that, we consider a coastal water with scintillation 
parameters of e = 10 _o m 2 /s 3 , w = —3 and \T = 4 x 10”' 
K 2 /s, and numerically evaluate Eq. ([3]) to find the scintillation 
index of a plane wave for different link ranges. Based on 
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Fig. 9. Comparison between the performance of a 30 m coastal water link 
with different configurations and a 25 m SISO link, both operating at R 5 = 2 
Gbps. 


Fig. 10. The exact BER of a SISO UWOC link for two separate scenarios 
(SI: a 25 m coastal water link with R 5 = 1 Gbps, and S2: a 10 m turbid 
harbor water link with Rt, = 200 Mbps) and four different models for the 
channel impulse response described within the text (Ml, M2, M3, and M4). 


our numerical results we find that for do = 25 m and 
30 m the log-amplitude variance a\ is 0.126 and 0.165, 
respectively. As it is obvious in Fig. 9, only 5 m (%20) 
increase on the communication range remarkably degrades the 
system performance, e.g., approximately 12 dB degradation 
is observed at the BER of 10~ 12 . But as it can be seen, 
increasing the number of independent links or equivalently 
mitigating fading deteriorations considerably improves the 
system performance. Since spatial diversity manifests itself 
as a reduction in fading variance 118 |. we can conclude that 
there exists a configuration with spatial diversity at link range 
of 30 m which performs similar to a SISO link at that range 
but with less fading variance, e.g., a\ = 0.01. Hence, in this 
figure the performance of a 30 m SISO link with <j\ = 0.01 
is also depicted for the sake of comparison and obviously 
it can yield better performance than a 25 m SISO link, 
especially for lower error rates. Therefore, one can achieve 
better performance even in longer link ranges by employing 
spatial diversity technique. We should emphasize that spatial 
diversity can provide more performance enhancement rather 
than those are presented in this paper; when the channel suffers 
from strong turbulence pO) . But since this paper is focused 
on weak oceanic turbulence, we only considered channels with 
a] < 1 (28). 

In order to better observe the effects of different factors 
on the channel impulse response and hence on the system 
performance, we consider four different channel models and 
evaluate the BER performance of a SISO UWOC link for 
each of the considered channel models. The first model (Ml) 
only considers turbulence effect and ignores absorption and 
scattering effects, i.e., hi(t) = a 2 6(t — d a /v), where v is 
the propagation speed of light through water. The second 
model (M2) only considers absorption and scattering using 
Beer’s law | 6 | and neglects turbulence effect, i.e., / 12 (f) = 
e~ cd °S(t — do/v). The third model (M3) obtains the channel 
fading-free impulse response using MC numerical simulations 
but ignores turbulence effect, i.e., / 13 (f) = ho(t). And the 
fourth model (M4), similar to the rest of the paper, takes all 
of the three impairing effects into account while considers 
absorption and scattering effects based on MC simulations and 
turbulence effect as a multiplicative fading coefficient, i.e.. 


/ 14 (f) = a 2 ho(t). We also consider two different scenarios: 
SI; a 25 m coastal water link with transmission rate Rb = 1 
Gbps, and S2; a 10 m turbid harbor water link with Rb = 200 
Mbps. All of the system parameters are the same as those are 
listed in Table I, and the log-amplitude variance for both of 
the scenarios is considered to be <j\ = 0.16. As it is shown 
in Fig. 10 based on the results of the accurate channel model, 
i.e., the fourth channel model M4, a 25 m coastal water link 
has approximately 18 dB better performance than a 10 m 
turbid harbor water link even for 5 times larger transmission 
rates. Furthermore, as expected, ignoring the turbulence effect 
through considering the channel impulse response as the 
third model M3, underestimates the channel impairments and 
significantly increases the slope of BER curves. On the other 
hand, considering the channel impulse response as the second 
model overestimates the channel attenuation and shifts the 
BER curves to the right side when compared with the third 
channel model results. This is because Beer’s law ignores 
those photons which may reach the receiver through multiple 
scattering. The overestimation of Beer’s law for turbid harbor 
water link is by far more than that of coastal water link 
due to the higher values of attenuation length defined as 
Tatn = cdo (for the above two scenarios T atn ,si = 9.95 and 
T~atn,S2 = 21.9) Q, (14). In other words, in turbid harbor 
waters the occurance of multiple scattering is very prevalent 
and many of the transmitted photons may reach the receiver 
after scattering many times, and ignoring multiple scattering in 
such links, through use of Beer’s law, results into an extreme 
overestimation of the channel loss. Finally, ignoring absorption 
and scattering effects and only considering turbulence effects 
as a multiplicative fading coefficient, through use of the first 
channel model Ml, significantly underestimates the channel 
degrading effects and considerably shifts the BER curves to the 
left side. This is mainly because in coastal and turbid harbor 
UWOC channels with weak turbulence conditions, absorption 
and scattering have much more impairing effects than turbu¬ 
lence, and ignoring their effects remarkably underestimates the 
channel impairments. 

As it has been observed, spatial diversity provides a sig¬ 
nificant performance improvement but under the assumption 
of independent links which is sometimes practically infeasible 
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O independent, SISO 

■ © — independent, 2x1 MISO 
O ' 1 independent, 3x1 MISO 
0 dependent, 2x1 MISO 

□ dependent, 3x1 MISO 

• dependent, 2x2 MIMO 


20 25 

Average transmitted power per bit [dBm] 


Fig. 11. Effect of spatial correlation on the performance of a 25 m coastal 
water link with ax = 0.4, Rb = 2 Gbps, p = 0.25 and 0.7. 


ED- Therefore, in practice received signals by different links 
may have correlation. In Fig. 11 we investigate the effect 
of spatial correlation on the performance of MIMO UWOC 
systems in a similar approach to 0- Given an M by M 
correlation matrix at the transmitter side as R r and also an N 
by N correlation matrix at the receiver side as Rr, the spatial 
correlation matrix of the MIMO channel can be obtained 
by the Kronecker product of the spatial correlation matrices 
of the transmitter and receiver, i.e., Rmimo = Rt®Rr, 
which is of size M x TV by M x N 0- Then, given 
X = [Xi, X 2 , ..., Xmxn\ as a vector of M x N Gaussian 
distributed RVs with mean zero and variance one, the vector 
of new RVs, X 1 = [X[, X' 2 ,X' MxN ], corresponding to 
the Gaussian distributed log-amplitude factors with correlation 
matrix of Rmimo. can be obtained as; 


X' = a x XC - £, 


(38) 


where S = [a \,..., cr\\ is a vector of M x TV elements of a\, 
and C is an M x TV by MxN upper triangular matrix obtained 
from the Cholesky decomposition of the correlation matrix 
Rmimo as Rmimo = C*C, where C* denotes the complex 
conjugate transpose of C. Finally, the new correlated lognormal 
RVs, corresponding to the fading coefficients of different links, 
can be obtained as a'^ = exp(2XC), and the average BER 
can be evaluated from these correlated fading coefficients by 
numerical calculation of multi-dimensional integrals. In the 


simulations of Fig. 11 we have considered a 25 m coastal 
water link with ax = 0.4 and lit, = 2 Gbps, and evaluated the 
upper bound BER of different configurations in several cases, 
i.e., independent links and correlated links with correlation 
values of b(lo) = p = 0.25 and 0.7. Comparing the results 
shows that the performance loss is much more severe for larger 
correlation values, e.g., 2 dB and 6 dB degradation can be 
observed in the performance of a 3 x 1 MISO link with the 
BER of 10 -10 and with p = 0.25 and 0.7, respectively. Also 
a 3 x 1 MISO link with p = 0.25 yields the same performance 
as an independent 2x1 MISO link, i.e., the diversity order is 
decreased by one. Therefore, the performance enhancement of 
spatial diversity in highly correlated weak turbulent channels 
may be insignificant. 


VI. Conclusion 

In this paper, we studied the performance of MIMO UWOC 
systems with OOK modulation and equal gain or optimal 
combiner. Our derivations were based on a general channel 
modeling which appropriately takes all of the channel im¬ 
pairments into account. In particular, we obtained the channel 
fading-free impulse response using MC numerical simulations 
and included turbulence effects as a multiplicative fading coef¬ 
ficient. Closed-form solutions for the system BER expressions 
obtained in the case of lognormal underwater fading channels, 
relying on Gauss-Hermite quadrature formula as well as 
approximation to the sum of lognormal random variables. We 
also applied photon-counting method to evaluate the system 
BER in the presence of shot noise. The excellent match 
between the results of analytical expressions and photon¬ 
counting method confirmed the validity of our assumptions 
in derivation of analytical expressions for the system BER. 
Furthermore, our numerical results indicated that EGC is more 
practically interesting, due to its lower complexity and its close 
performance to optimal combiner. In addition to evaluating 
the exact BER, also the upper bound on the system BER 
has been evaluated and excellent tightness between the exact 
and upper bound BER curves has been observed. Moreover, 
the good match between the results of numerical simulations 
and analytical expressions verified the accuracy of our derived 
analytical expressions for the system BER. Our numerical 
results showed that spatial diversity manifests its effect as 
a reduction in fading variance and hence can significantly 
improve the system performance and increase the viable 
communication range. In particular, a 3x 1 MISO transmission 
in a 25 m coastal water link with log-amplitude variance of 
0.16 can introduce 8 dB performance improvement at the BER 
of 10 9 . We also observed that spatial correlation can impose 
a severe loss on the performance of MIMO UWOC systems. 
Specifically, correlation value of p = 0.25 between the links 
of a 3 x 1 MISO UWOC system with ax = 0.4 decreases the 
order of diversity by one. Finally, we should emphasize that 
although all of the numerical results of this paper are based on 
lognormal distribution, many of our derivations can be used 
for any other fading statistical distribution. 


Appendix A 

Decision Rule for MIMO UWOC System with OC 
In this appendix, we obtain the decision rule for BER 
evaluation of MIMO UWOC system with OC. Based on the 
ML detection rule and assuming the availability of perfect 
CSI, the symbol-by-symbol receiver which does not have any 
knowledge to {6; c } fc " = _ / , adopts the following metric for 

optimum combining 0; 

Pr(r|6 0 = l,a)^Pr(f|6 0 = 0, a), (39) 

0 

where f = (rr . 72 ,.... rx) is the vector of different branches’ 
integrated received current and a = (an, ai 2 ,..., aMiv) is the 
fading coefficients vector. The conditional probabilities for the 
“ON” and “OFF” states are respectively given as; 


Pr(r|6 0 = !,«) = 


1 


-1 


(27r4J JV/: 


exp 


N 

Z<Jt » 3 = 1 


M 




i =1 


(40) 
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Pr (r|&o = 0, a) = 


(2 ™l) N /‘ 


exp 


N 

E 5 


-1 

2(7™ J 

T,, i=i 


■ (41) 


Replacing ( |40| and ( |4T] > in ( |39| ) and dropping the common 
terms out, the decision rule simplifies to; 

2 

I ■ (42) 


N M N 


M 


EEtf 




j=i 


3= 1 \<=1 


Appendix B 

(M x 7V)-Dimensional Series of Gauss-Hermite 
Quadrature Formula 

In this appendix, we show how (M x TV)-dimensional 
averaging integrals over fading coefficients can effectively 
be calculated using Gauss-Hermite quadrature formula [37, 
Eq. (25.4.46)]. More specifically, we calculate = 

Is bfc /< 3 (cE)da with (M x ^-dimensional series 

where ^ is defined in ( | 1 4[ > and ( | 1 5[ ) (e.g., for OC) for 

bo = 1 and bo = 0, respectively. Based on ( [T2| ) and [37, Eq. 
(25.4.46)], for any function g(afj) averaging over lognormal 
distributed fading coefficient a. l3 can be calculated with a finite 
series as; 

/*oo 

/ dfaij^fanj {p^ij)dOLij ~ 

JO 

■4 E = ex P ( 2 4oV 2CT i- + 2 f^.,) )■ 

V 5ij=l 

(43) 

Further, the validity of Eq. ( |44| ), shown at the top of 
the next page, can be verified by induction for any func¬ 
tion of M x N lognormal distributed fading coefficients 
g{a\ 3 , c* 2 i,cx 2 MN ). Therefore, the (M x tV)-dimensional 
integral of can be calculated as Eq. 

m 


Appendix C 

MGF of the Receiver Output in SISO Scheme 

In this appendix, we calculate the receiver output MGF in 
SISO scheme. Based on (|27j), conditioned on {b k } k f_ L and a, 


,(M 

x SISO 


is the sum of two independent RVs. Therefore its MGF 
T (6 o) (s) is the product of their MGFs, i.e., T ,b 0 ) (s) = 

^SISO ^SISO 

T (b 0 ) (s) x & Vth (s). We first obtain the conditional MGF of 

^SISO 

t/gjgQ conditioned on {b k } k= _ L and a. Then averaging over 
{bk}kf-L resu l ts the MGF of t/gjgo conditioned on a 136 


Eq. ( |46l ). shown at the top of the next page page. Note that 
b k s are independent Bernoulli RVs with identical probability, 
i.e., Pb k {bk) = \S{b k ) + \5(bk — 1). Therefore, the latter 
expectation in (|46| simplifies to; 


E, 




-1 

[ [] exp (b k a 2 m^^ (e s - 1) 


k——L 
-1 


= n e * 

k=—L 
-1 

= n (V 2 ) 


exp 


[bkd 


2 m (/A) 


(e s - 1) 


1 + exp ( a 2 m < - I ’ k ' > (e s — 1 


(47) 


fc=-L 


Finally, inserting ( |47] > in ( |46| ) and then multiplying the result 
by = exp(s 2 (7j /l /2) yields the output MGF as in Eq. 
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Appendix D 

MGF of the Receiver Output in MIMO Scheme 


In this appendix, we calculate the MGF of the receiver 
output in MIMO scheme. Since EGC is used, the combined 

Conditioned 


output of the receiver is i%imo = Syl 
on fading coefficients vector a, received signals from different 
branches are independent and hence MGF of their sum is 

the product of each branch’s MGF, i.e., _ (s) = 

“mimo 1“ 

we first need to obtain 


nf=i* 


uf 0) Haijl^i. 


(s). Therefore, 


the conditional MGF of each branch. By pursuing similar 


steps as in Appendix C, ’T 


Oo) 


as Eq. 




(s) can be calculated 


as Lij = L 


shown at the top of the next page, where 
= (ribj +n ( ij)T} > . Supposing the same channel memory 
iax for all links, performing the latter expectation 
in ( |48j ) results the jth receiver output MGF as; 

/ 2 9 T m 

( S )=ex pf-^| 




(<>o) 


I kj}“ 


(bd) 


^bpa 


2 (s) 

iMJ 


(e s -l) 


II (V 2 ) 


k——L n 


M 

1 +n exp ( 

2=1 


2 (/,fc) 


(e s -l 


(49) 


Note that if some of the links have smaller memory, we can 
add adequate zero components to each of them in order to 
make all the links with an equal memory length. In other 


words, if La < L max , then {m\ ■ } k= l= 0. Finally, 


-‘V "" ^max, \"H,j ffc=-L max 

MGF of the receiver output in MIMO scheme can be obtained 
as in Eq. 
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